/*
 * raytracer.cu - raytracing gpu code
 * Wingfai Chan, Daniel Kotowski
 * CIS 4450
 *
 */


#include <vector>
#include <vector_types.h>
#include <vector_functions.h>
#include <cutil_math.h>

texture<float4,1,cudaReadModeElementType> triangles_texture;
texture<float4,1,cudaReadModeElementType> bounding_texture;

/* 
 * Convert floating point rgb color to 8-bit integer
 */
__device__ int rgbToInt(float r, float g, float b)
{
	r = clamp(r, 0.0f, 255.0f);
	g = clamp(g, 0.0f, 255.0f);
	b = clamp(b, 0.0f, 255.0f);
	return (int(r)<<16) | (int(g)<<8) | int(b); // notice switch red and blue to counter the GL_BGRA
} 

/*
 * Calculates the distance between two points in euclidean coordinates
 */
extern "C" __device__ __host__ float dist(float3 point1, float3 point2)
{
	return (sqrt(powf(point1.x-point2.x, 2) + powf(point1.y-point2.y, 2) + powf(point1.z-point2.z, 2)));
}

/*
 * Determine if a ray intersects a triangle and returns the distance; based on the implementation found in this paper
 * http://www.cs.virginia.edu/~gfx/Courses/2003/ImageSynthesis/papers/Acceleration/Fast%20MinimumStorage%20RayTriangle%20Intersection.pdf
 */
__device__ float rayTriangleIntersect(float3 ray_origin, float3 ray_dir, float3 vertex0, float3 edge1, float3 edge2, float &u, float &v)
{
	float3 tvec = ray_origin - vertex0;  
	float3 pvec = cross(ray_dir, edge2);  
	float  det  = dot(edge1, pvec);  
	float inv_det = __fdividef(1.0f, det);  

	u = dot(tvec, pvec) * inv_det;  

	if (u < 0.0f || u > 1.0f)  
		return -1.0f;  

	float3 qvec = cross(tvec, edge1);  

	v = dot(ray_dir, qvec) * inv_det;  

	if (v < 0.0f || (u + v) > 1.0f)  
		return -1.0f;  

	return (dot(edge2, qvec) * inv_det);
}

/*
 * Determines if the sphere intersects the ray
 */
__device__ float raySphereIntersect(float3 sphere_centre, float sphere_radius, float3 ray_origin, float3 ray_dir)
{
	
	float3 dst = ray_origin - sphere_centre;

	float A = dot(ray_dir, ray_dir);
	float B = 2*dot(dst, ray_dir);
	float C = dot(dst, dst) - sphere_radius*sphere_radius;

	//check the determinant to see if there is an intersection
	if (B*B - 4*A*C > 0)
		return 1;
	else
		return 0;
	

/*
	float3 dst = ray_origin - sphere_centre;

	float B = dot(dst, ray_dir);
	float C = dot(dst, dst) - sphere_radius;

	//check the determinant to see if there is an intersection
	if (B*B - C > 0)
		return 1;
	else
		return 0;*/
}

/* Calculates if the object is in shadows*/
__device__ int inShadow(float3 origin, float3 ray, int num_bounding) 
{
	int i, j;
	float4 tmp;
	float3 vertex, edgea, edgeb;
	float t_intersect, u, v;

	
	//go through the list of bounding volumes
	for(j=0; j<num_bounding; j+=2)
	{
		//if we intersect one
		tmp = tex1Dfetch(bounding_texture, j);
		if (raySphereIntersect(make_float3(tmp), tmp.w, origin, ray))
		{
			//get the start and end index of the triangles and test for intersection
			tmp = tex1Dfetch(bounding_texture, j+1);
			for (i=(int)tmp.x; i<=(int)tmp.y; i+=7)
			{
				vertex = make_float3( tex1Dfetch(triangles_texture,i));
				edgea  = make_float3( tex1Dfetch(triangles_texture,i+1));
				edgeb  = make_float3( tex1Dfetch(triangles_texture,i+2));
				
				t_intersect = rayTriangleIntersect(origin, ray, vertex, edgea, edgeb, u, v);
				
				if (t_intersect > 0.001)
					return (1);
			}
		}
	}

	return(0);
}

__device__ float3 hitObject(float3 origin, float3 ray, int num_bounding,float3 light_pos,float3 light_colour)
{
	int i, j, vnum,depth;
	//Fetched Variables
	float3 vertex;
	float3 edgea;
	float3 edgeb;
	float4 tmp;
	float3 normal0;
	float3 normal1;
	float3 normal2;

	float t_intersect,reflct, u, v, closest_u, closest_v;
	float intersect;
	float3 colour;
	float3 base_colour;
	float3 hit_position;	
	float3 vec_to_light;
	float3 point_normal;
	float3 reflection;
	float3 specular;
	float DOT;

	colour.x = 0;
	colour.y = 0;
	colour.z = 0;
	
	//loop for calculating reflections since CUDA doesn't allow recursion
	for (depth=0; depth<4; depth++)
	{
		intersect = UINT_MAX;
		/*Calculate Closest triangle*/
		for(j=0; j<num_bounding; j+=2)
		{
			tmp = tex1Dfetch(bounding_texture, j);
			if (raySphereIntersect(make_float3(tmp), tmp.w, origin, ray))
			{
				//return (make_float3(255, 0, 0));
				tmp = tex1Dfetch(bounding_texture, j+1);
				for (i=(unsigned int)tmp.x; i<=(unsigned int)tmp.y; i+=7)
				{
					vertex = make_float3( tex1Dfetch(triangles_texture,i));
					edgea  = make_float3( tex1Dfetch(triangles_texture,i+1));
					edgeb  = make_float3( tex1Dfetch(triangles_texture,i+2));
					
					t_intersect =	rayTriangleIntersect(origin, ray, vertex, edgea, edgeb, u, v);
					
					if(t_intersect >= 0.001 && t_intersect < intersect)
					{
							intersect = t_intersect;
							vnum = i;
							closest_u = u;
							closest_v = v;
					}
				}
			}
		}
			
		if(intersect != UINT_MAX)
		{
			
			vertex = make_float3(tex1Dfetch(triangles_texture, vnum));
			edgea  = make_float3( tex1Dfetch(triangles_texture,vnum+1));
			edgeb  = make_float3( tex1Dfetch(triangles_texture,vnum+2));
			normal0 = make_float3(tex1Dfetch(triangles_texture,vnum+4));
			normal1 = make_float3(tex1Dfetch(triangles_texture,vnum+5));
			normal2 = make_float3(tex1Dfetch(triangles_texture,vnum+6));

			tmp =  tex1Dfetch(triangles_texture,vnum+3);		
			reflct = tmp.w;
			base_colour = make_float3(tmp);
			
			//calculate specular and diffuse lighting
			hit_position = intersect*ray+origin;
			vec_to_light = normalize(light_pos - hit_position);	

			//point normal interpolated using the barycentric coordinates from the intersection
			point_normal = normalize(normal0 + closest_u * (normal1 - normal0) + closest_v * (normal2 - normal0));
			//point_normal = plane_normal;
			reflection = normalize(reflect(vec_to_light,point_normal));
			DOT = max(dot(ray,reflection),0.0);

			specular =__powf(DOT,25)*make_float3(255,255,255);
			
			colour += base_colour*0.2 + base_colour*max(dot(point_normal,vec_to_light),0.0)+specular*reflct;
			
			//calculate shadows
			if(inShadow(hit_position,vec_to_light, num_bounding)==1)
			{
				colour *=.6;
			}
			
			//calculate reflections
			if(reflct==1.0)
			{
				origin = hit_position;
				ray  = normalize(reflect(ray, point_normal));
				colour = make_float3(0, 0, 0);
			}
			else
				break;
		}
		else
		{
			colour += make_float3(135, 206, 250);
			break;
		}
	}

	if (depth>0)
	{
		colour /= depth;
	}
	/*

		return(colour);
	}*/

	/* For each triangle
		  if ray intersects 
			  if Closest intersection 
		if(intersected)
			return interseciton point and normal
	*/

	return(colour);

}

void __global__ rayTrace(unsigned int *out_data, int width, int height, int num_volumes, float3 camera_pos, float3 a, float3 b, float3 c, float3 light_pos, float3 light_colour)
{
	int x = blockIdx.x * blockDim.x + threadIdx.x;
	int y = blockIdx.y * blockDim.y + threadIdx.y;
	
	float xf = (x-0.5)/((float)width);
	float yf = (y-0.5)/((float)height);

	/* Calculate Ray for specific pixel */
	float3 ray;
	float3 colour;

	float3 t1 = c+(a*xf);
	float3 t2 = b*yf;
	float3 image_pos = t1 + t2;
	ray = normalize(image_pos-camera_pos);
	//ray = image_pos-camera_pos;
	
	
	if (raySphereIntersect(light_pos, 2.0, image_pos, ray))
		colour = make_float3(255,255,255);
	else
		colour = hitObject(image_pos,ray,num_volumes,light_pos,light_colour);
	out_data[y*width + x] = rgbToInt(colour.x, colour.y,colour.z);
}

extern "C" void launchRayTrace(unsigned int *out_data, int width, int height, int num_volumes, float3 camera_pos, float3 a, float3 b, float3 c, float3 light_pos, float3 light_colour)
{
	dim3 blocks(width/8, height/8, 1);
	dim3 per_block(8, 8, 1);
	rayTrace <<<blocks, per_block>>>(out_data, width, height, num_volumes, camera_pos, a, b, c, light_pos, light_colour);
}

/*
 * This function binds the triangles on the device as a 1D texture reference
 */
extern "C" void bindTriangleTexture(float *triangles, size_t size)
{
	//bind the texture memory
	triangles_texture.normalized = false;
	triangles_texture.filterMode = cudaFilterModePoint;
	triangles_texture.addressMode[0] = cudaAddressModeWrap;

	cudaChannelFormatDesc channel_desc = cudaCreateChannelDesc<float4>();
	cudaBindTexture(0, triangles_texture, triangles, channel_desc, size);
}

/*
 * This function binds the bounding volumes on the device as a 1D texture reference
 */
extern "C" void bindBoundingTexture(float *bounding_volumes, size_t size)
{
	//bind the texture memory
	bounding_texture.normalized = false;
	bounding_texture.filterMode = cudaFilterModePoint;
	bounding_texture.addressMode[0] = cudaAddressModeWrap;

	cudaChannelFormatDesc channel_desc = cudaCreateChannelDesc<float4>();
	cudaBindTexture(0, bounding_texture, bounding_volumes, channel_desc, size);
}